i/1 


I*  AD-R132  885 
UNCLASSIFIED 


NUMERICAL  ANALYSIS  OF  BOUNDARY  VALUE  PROBLEM  OF 
ELLIPTIC  TYPE  BY  MEANS  PE.  (U)  WISCONSIN  UNIV-MADISON 
MATHEMATICS  RESEARCH  CENTER  T  HANADA  ET  AL.  JUL  83 
MRC-TSR-2541  DAAG29-80-C-0041  F/G  12/1  NL 


MRC  Technical  Summary  Report  #  2541 


Mathematics  Research  Center 
University  of  Wisconsin-Madison 
610  Walnut  Street 
Madison,  Wisconsin  53706 


July  1983 


(Received  May  25, 


1983) 


m  FILE  COPY 


Sponsored  by 

U.  S.  Army  Research  Office 
P.  O.  Box  12211 
Research  Triangle  Park 
North  Carolina  27709 


Approved  for  public  release 
Distribution  unlimited 


UNIVERSITY  OP  WISCONSIN  -  MADISON 
MATHEMATICS  RESEARCH  CENTER 


NUMERICAL  ANALYSIS  OP  BOUNDARY  VALUE  PROBLEM  OF 
ELLIPTIC  TYPE  BY  MEANS  PENALTY  AND  THE 
FINITE  DIFFERENCE  AND  ITS  APPLICATION  TO  FREE  BOUNDARY  PROBLEM 

*  •*  *** 

T.  Hanada  ,  H.  Kawarada  and  O.  Pironneau 

Technical  Summary  Report  #2541 
\  July  1983 

_  ABSTRACT 

\  T/<t 

— ^w^  study  a  numerical  method  for  solving  free  boundary  problems  of 
elliptic  type.  Usually  these  problems  are  prescribed  with  two  boundary 
conditions  on  the  free  boundary.  One  of  them  is  the  Dirichlet  condition  and 
the  other  is  the  Neumann  condition.  Our'  method  is  to  transform  the  original 
problem  to  an  optimization  problem.  The  state  equation  is  approximated  by  an 
equation  with  a  penalty  terra  in  which  the  Dirichlet  condition  on  the  free 
boundary  is  approximately  satisfied.  The  outward  normal  derivative  included 
in  the  Neumann  condition  through  the  free  boundary  is  calculated  by  using  the 
asymptotic  behavior  of  the  solution  of  the  penalized  state  equation,  -ncre  we 
— present-  a  method  to  solve  this  penalized  optimization  problem.  Also  the  error 
estimate  of  the  discretized  state  equation  by  the  finite  difference  method  is 
given. 
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SIGNIFICANCE  AND  EXPLANATION 


Free  boundary  problems  for  elliptic  partial  differential  equations  arise 
in  many  applications,  e.g.,  phase  transition  problems  including  slag  flow  in 
the  hearth,  a  wave  or  jet  problem,  the  equilibrium  of  plasma,  optimal  shape 
design  and  others.  Hence  there  is  interest  in  the  development  of  efficient 
and  accurate  numerical  methods  for  the  solution  of  these  problems. 

In  this  report  we  develop  a  method  for  the  numerical  treatment  of  these 
problems.  In  this  method  we  replace  the  original  problem  by  an  optimization 
problem  in  which  the  state  equation  contains  a  "penalty"  term.  This  is 
particularly  useful  because  the  boundary  conditions  on  the  free  boundary  are 
satisfied  approximately;  this  avoids  difficult  calculations.  Another 
objective  of  this  report  is  to  give  the  error  estimate  of  the  solution  of  the 
discretized  state  equation  by  the  finite  difference  method.  Our  proof  uses 
some  of  the  techniques  of  the  finite  element  method;  it  requires  a  new 
estimate  of  the  solution  of  the  penalized  state  equation. 
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1.  Introduction 

This  paper  is  concerned  with  a  numerical  method  for  solving 

2 

the  problem  of  Dirichlet  in  an  arbitrary  domain  of  R  using 
the  finite  differences.  Consider,  for  example,  the  problem 
of  finding  <J>°  such  that 


(1.1) 


-  A  i|>°  +  Xf  i|>°  =  f 


in  ft0. 


*°l  o 


=  0  . 


Let  Q  be  a  rectangle  which  includes  I)  and  let  x  be  the 
characteristic  function  of  then  (1.1)  may  be  approximated 


(1.2) 


-a»i>g  +  xn<i>G  +  -  (i-x)f 


in  £2  , 


<|>e|3£2  =  0  . 


This  penalization  method  was  studied  in  Kawarada  (1977) ,  (1979) 
where  it  was  shown  that 


(1.3) 


/F  I 
a  SSI0 


+  0(e) 


Department  of  Information  Mathematics,  The  University  of  Electro 
Communications . 

*  *  , 

Department  of  Applied  Physics,  University  of  Tokyo. 

***  _  , 

Department  of  Mathematics,  University  of  Pans-Nord. 


When  (1.2)  is  approximated  by  using,  say  the  5  points  for¬ 


mula  for  the  Laplace  operator,  then  y  must  be  approximated  in 
such  a  way  as  to  yield  a  computable  linear  system  for  the  values 
of  \|ie  at  the  grid  points.  The  basic  difficulty  with  such  method 
is  that  it  is  not  easy  to  obtain  a  good  approximation  of  3$C/3n 
on  3ft®  from  the  discrete  approximation  of  $c. 

In  this  paper  we  shall  show  that  x  can  be  approximated 
(see  (2.14)- (2.15) )  so  as  to  yield  (Theorem  1  and  4). 

d.4)  ||$  5  -  $e  !|  «  o(h2/3)  <*> 

n  i,n 


(1.5) 


(-  -  $5) 

/l  h 


l/2,3ft° 


0(h1/3) 


» 


where  h  is  the  mesh  size,  e=h2/3,  the  norms  are  Sobolev 

norms,  and  <f>E  is  the  linear  interpolate  from  the  values  of  the 

approximation  of  computed  at  the  grid  points. 

The  proof  uses  some  of  the  techniques  of  the  finite  element 

method  (see  Ciarlet  (1978)  for  example);  it  requires  the 

estimates  for  ||$e-$°||_?  <*>  and  ll*e|L2  »  which  were  not  known 

11  w2'  (n0)  w2'  (n-n0) 

before;  thus  this  paper  contributes  also  to  the  theory  of  asymp¬ 
totic  expansion  in  this  respect. 

The  method  of  penalization  used  here  (see  (1.2))  belongs 
to  the  family  of  artificial  domains,  capacitance  matrices ...  for 
which  the  scientific  literature  abonds  (see  for  example  Lions- 
Marchuck  (1979),  Proskurowski  and  Widlund  (1966),  (1980)  and  the  biblio¬ 
graphy  there  in) ;  it  is  quite  possible  that  the  techniques  used  here  apply 
also  to  these  other  methods.  Finally  our  practical  purpose  is  to 

^  ^  II vll  _  indicates  the  norm  of  v  in  Hm . 


■m 

j 


obtain  the  solution  of  a  free  boundary  problem  transformed  by 
least  squares  into  an  optimum  design  problem;  the  method  goes 
as  follows: 


Let  be  bounded  by  a  given  boundary  Tq 
y  (cf.  Fig.l). 

We  are  looking  for  $  and  Y  such  that 


and  a  free  boundary 


(1.6)  -  A  +  +  XQ*  -  f  in  ftY 

ci.7)  *|r  -  o  ,  r  *  «  r0  u  y  (r0  n  y  -  ❖> 


(1.8) 


2 

where  f  and  g  are  given  in  R  and  is  a  positive  constant, 
n  is  outward  normal  to 


Y 


Figure  1 

We  can  now  transform  this  problem  into  a  problem  of  optimum 
design*  which  is  an  optimal  control  problem  of  a  distributed 
parameter  system  where  the  control  is  a  part  of  the  boundary. 


For  a  given  y,  we  define  the  state  equation  by 


(1.9)  ♦  X0*  «  f  in  fty 

(1.10)  *|r  -  0 

which  defines  4»  *  4>(x;Y)  •  Then  any  solution  y  of  the  problem 
(1.6)-(1.8)  is  also  solution  of 

(1.11)  min  (E(y)  -  f  .  g|2  ds}  . 

res  *  y  3n 

where  S  is  the  set  of  admissible  boundaries  y  such  that  3fl^  =  r^u  Y 
and  the  state  equation  has  solution  in  Problems  of  the 

type  of  (1.11)  was  successfully  studied  by  Murat-Simon  (1977), 
Chesnais  (1975)  ,  Pironneau  (.1976)  and  Deryieux  (1981) . 

The  characteristic  of  our  method  to  solve  (1.11)  is  sum¬ 
marized  as  follows: 

1*  Let  Q  be  a  bounded  domain  in  R*  such  that  end 

(I^nr^#)  for  VyfcS  (cf.  Fig.2).  We  penalize  the 
state  equation  as  follows:  For  e>0, 

(1.12)  (-A  +  X0He  +  |  Xy4»E  -  d-Xy)f  in  « 


(1.13) 


where  y  is  the  characteristic  function  of  (l  -  Q  .  If  e  is  small 

Ay  y 

enough,  it  is  shown  that  the  solution  4<e  of  (1.12)  and  (1.13) 
has  the  following  asymptotic  behavior  on  y: 


(1.14) 


+  0(e) 


in  a  suitable  topology 


* 


m 


■  .*> 
.* ; 


i 


n 
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where  is  the  solution  of  the  state  equation.  This  property 
will  be  proved  in  the  next  section.  From  (1.14)  we  observe 
that  the  Dirichlet  boundary  condition  (1.13)  is  approximately 

is  approximated 


satisfied  in  the  error  of  0(/e)  and  also  ~ 
.  on 

..  e « 


by  -  xK 

method. 


in  the  same  order,  which  constitutes  the  key  of  our 


2*  By  using  (1.14),  we  approximate  (1.11)  in  the  following 
royt 


1.15) 


min  (E£(y) 


f  +  g|2  ds) 


3* 


After  the  penalized  state  equation  (1.12)  and  (1.13) 


is  discretized  by  the  method  of  finite  difference,  (1.15)  will 
be  solved  by  the  numerical  technique  for  solving  problems  of 
optimum  design  (Pironneau  (1983)). 

Finally  we  would  like  to  mention  about  the  discretization 
of  the  penalized  state  equation  in  3s.  In  order  to  solve  the 
optimazation  problem  like  (1.11)  or  (1.15),  usually  the  gradient 
method  is  applied  together  with  the  finite  element  method;  in 
which  case,  the  mesh  is  moved  at  each  iterative  approximation 
of  the  free  boundary.  We  have  frequently  experienced  that  these 
moving  mesh  causes  serious  computational  problems. 

The  plan  is  following: 

1.  Introduction 

2.  The  Penalty  Method  (including  main  theorems) 

3.  Proof  of  Theorem  2 

4.  Proof  of  Theorem  3 

5.  Proof  of  Theorem  4 

6.  The  algorithm  for  solving  free  boundary  problem 
2.  The  Penalty  Method 

The  aim  of  this  section  is  to  give  the  proof  of  (1.5), 

2.1  Let  us  consider  the  following  problem;  given  fig  with  suffi¬ 
ciently  smooth  boundary  r,  (JTq<=  ft  =  (0. 1) x(o.  1 ) )  find  such 
that 

(2.1)  “  A  =  f  in 

(2.2)  *0L  -  0  ,  r  -  30*  . 


0 

.  1 

:■! 

* 

_  « 


-  1 


Then  the  solution  of  this  problem  satisfies 

(2.3)  *°  €  Hra+1(n0)  D  Hj(nQ) 

«•« 

for  any  ft  (m>0). 

2.2  We  shall  penalize  the  problem  (2.1)  and  (2.2)  as  follows; 
For  >  o, 

(2.5)  (-A  +  X0)*e  +  ^X<I»C  -  (l-X)f  in  « 

(2.6)  *€Ua"° 

where  x  is  the  characteristic  function  of  n^*n-TTg.  It  is 
checked  that  the  solution  <pC  satisfies 

(2.7)  *e  e  H2(Q)  H  Hj(fi) 

(2.6)  *►  (1-X)^^  strongly  in  H*(Q)  as  e  +  0. 

2.3  Let  us  approximate  (|>e  by  the  solution  of  the  discrete 
problem  corresponding  to  (2.5)  and  (2.6): 

<<-*„  +  +  !«vE>ij  - ((1  -  xh'f>ij  • 


(2.9) 


ft  is  the  set  of  lattice  points  which  divide  ft  into  N  *  N  equidis¬ 
tant  meshes  3ft  consists  of  the  lattice  points  laid 

on  3ft.  is  the  five  points  formula  fcr  A.  In  order  to  define 

the  modified  characteristic  function  xh»  we  need  some  notations 
and  definitions. 

(2.11)  *  ^Wh  ^  is  affine  on  } 

Tj «*  triangle  obtained  by  dividing  each  rectangle  of  the  mesh  by 
its  first  diagonal  (cf.  Fig. 3). 


Figure  3 


Define  (x)  by 

^  H^(ft)  and 

(lh  ,  mh)**^  for  any  (1  ,  m)  . 

Further  define 

(2.12)  o*j  =  (x  |  wlj(x)  >  0} 

n 


-8- 


(2.13) 


2^  -  n  no  *  ♦  and  °h^  A  Qi  * 


Then 

(2.14) 


Uj>*fch 


h 


defines  a  narrow  strip  along  r.  With  the  notations  above. 


*h- 

Xh<x)  is 

defined 

so  as  to  satisfy  the 

following  conditions 

i) 

Xh(x)  -  1 

in 

fllh  “  ni  "  ni  n  Dh? 

ii) 

Xh(x)  =  0 

in 

fl0h  “  n0“n0  n  Dh' 

XL) 

0  <  Xh(x) 

<  1 

in  Dh 

iv) 

H  RS  ||  - 

^(s)^ 

(x) dr1  (see  Pig.  4,  ||RS||  implies  the 

distance  between  R  and  S) . 

The  meaning  of  iv)  is  clarified  in  the  proof  of  Theorem  4  stated 
later.  On  the  other  hand,  it  plays  an  important  role  in  solving 
(1.15). 


Figure  4 


As  two  examples  of  x^/  we  give 


(2.15) 


Xh(x) 


^(x) 


in  Dv, 


(2.16) 


Xh(x) 


p+*-2(x) 


in  D. 


where  Jl1  (x)  and  i2  (x)  are  constant  on  r1(S)  and  equal  to 
||RS||  and  ||TS||  respectively,  and  p  =  ||x-S||. 

The  linear  interpolation  of  ((i»j)6&)  is  represented  by 


(2.17) 


♦h<*>  -  l 


2.4  Now  we  shall  state  the  following  approximate  theorem  for 
■^_|r  which  plays  an  essential  role  in  our  method: 


Theorem  1 
small  enough. 
Then 


Suppose  f  6  Hm(nQ)  (m>5)  and  e  **  h2^3  and  let  h  be 


(2.18) 


^  +J  ♦£  Hi.  r  *  °(hJ>  ' 


In  order  to  prove  Theorem  1,  we  have  to  prepare  the  follow¬ 
ing  three  theorems: 

Theorem  2  Suppose  ffcH^fln)  (m^O)  and  let  e  be  small  enough. 


2 


e 


Theorem  3 


Suppose  f  e  H?l(nQ)  (m  £  5)  and  let  e  be  small  enough 


Then 

(2.20) 

ii  *‘nt 

m 

(Q)  0 

(2.21) 

ii  yt 

I  < 

C1(l+/e) 

It  ?xpxj 

c<V 

(2.22) 

II  a2*' 

f2 

||  3xpxJ 

cor.) 

/e 

(i  =  0,  1,  2)  is  a  positive  constant  independent  of  e. 

The  regularity  properties  in  Theorem  3  are  needed  to  prove 
the  following  theorem. 


Theorem  4  Suppose  f€Hra(flQ)  (in  £5)*  Then 


(2.23) 


(»  0(hJ)  when  e«h  j 


2.5  Here  we  give  the  proof  of  Theorem  1.  From  Theorems  2  and 
4,  we  have 


(2.24) 


(2.25) 


3i|>° 

"In" 


1 

7 


r 


0(/e) . 


Combining  (2.24)  and  (2.25) , 


win 


(2.26) 


♦h 

✓c 


(2.27) 


/c 


3i|>0 

^n" 

1  r 
7'  r 

h2'3. 

then 

»♦•  ii 

’  In 

Hi.r 

0(/c+|+  (|)2). 


1 

0(h7) 


3.  Proof  of  Theorem  2 

3.1  We  introduce  some  operators  defined  between  traces  on  T. 

(i)  Define  the  mapping 

\  3*al  *  7 

(3.1)  Tf  *  H  (r)  9  a  ~ -§ir|r e  H  tr)  : 

is  the  solution  of  the  problem: 

(3.2)  -A*  +  A0*  «  f  in  «0 

(3.3)  *lr  «  a 

where  f  6  h”"”1  (Dq)  . 

(ii)  Define  the  mapping 

-  1  1 

(3.4)  Re  :  H  7(D  B  b  «•  *5 1  6  H2(T)  : 

D|r 

is  the  solution  of  (3.2)  with  f  5  0  and  the  boundary  condition 

t  D 


-12- 


(lii)  Define  the  mapping 


(3.6) 


Sc  :  I?(D  9  a 


-  T 

!F|reH  (F) 


where  is  the  solution  of  the  problem: 

O 


(3.7) 


(3.8) 


(3.9) 


-eA\p  +  ip  *  0  in  ft. 


*lr  -  » 


♦laa  *  0  * 


We  denote  T^f  Se  and  Re  by  the  restriction  of  T*,  Se  and  R 
1  r  in  ra  r 

m+  -y 

to  H  (D .  But,  we  abbreviate  the  suffix  m  hereafter. 


Lemma  3.1  Let  a,  b  be  arbitrary  in  h""1'2  tr> .  Then 


(3.10)  Tf(a)  -  Tf(b)  »  TQ(a-b) 


where  TQ  *  Tf . 

Proof  Let  <ji„  (s=a,b)  be  the  solution  of  (3.2)  and  the  bound- 
* s 


ary  condition 


(3.11) 


*lr  *s  . 


Put  ¥  *  7  satisfies 


(3.12) 


-A7  +  -  0  in  ft„  , 


(3.13) 


»|r  -  a  -  b 


From  (3.14)  and  (3.15)  follows  (3.10).  | 

Here  it  should  be  noted  that  TQ  is  linear  and  Tf  is  non-linear. 


Lemma  3.2 

Tf,  (Rc)_1  and  Seare  homeomorphic  from  Hn“1/,2(D  to  Hm+1/^2(r)(m  >  Q) , 
Proof 

1*  Tf  is  injective  from  Hm+1/2 (D  into  Hm"1/2(D.  Indeed,  let 
a,b*Hm+1/2(r)  (a*b).  Suppose  Tf  (a)  -  Tf  (b) .  Then,  by  (3.10) 

0  -  Tf  (a)  -  Tf  (b)  *  TQ(a-b)  *  0 

because  of  the  strong  maximum  principle  under  the  assumption 
Xq  >  0.  This  is  a  contradiction. 

2*  Tf  is  surjective  from  Hm*1/2(r)  onto  Hm"1/2(D  if  XQ  >  0. 

Zn  fact,  choose  any  b  €  H*0”^^2  (D  .  Then  the  following  problem: 

(3.16)  -A*  +  X0i|>  «  f  in  nQ 


(3.17) 


has  a  unique  solution  Hm+^‘ (Qq)  if  XQ  >  0,  which  satisfies 


-14- 


3*  It  is  checked  that  Tf  and  (Tf)  are  continuous  between 
Hm+1/^2(D  and  Hm"^2(D  (see  Agraon,  Douglis  and  Nirenberg  (1959)). 

4*  Summing  up  1*,  2*  and  2*,  we  see  that  is  a  homeomorphism 
from  Hm+1/2(D  onto  Hm"1/2(D. 

The  repeated  use  of  the  above  arguments  proves  that  Re  and  S 
are  homeomorphic  between  Hm+1/2(n  and  H““1/2(r). 

■ 

3.3  Here  we  give  the  estimates  of  the  norm  of  Re  and  Se, 
which  are  crucial  for  the  proof  of  Theorem  2. 

Lemma  3.3  Let  e  be  small  enough  and  m^O.  Then 

(3.19)  ||  R€  (a)  ||  ,  0(1)  ||  a  ||  ,  ,  for  Va  e  ^1/2  (r 

s*  j #  r  b« ^ ,  r 

(3.20)  ||  Re  (a)  ||  ,  .  0(  i)  ||a  |j  .  ,  for  Va  e  H m"1/2  ( 

b  + j »  r  /e  b*  j  i  r 


(3.21)  ||  Rc  (a)  -a  ||  ,  «  o(/c)  ||  a  ||  , 

a- j »  r  b+  j  *  r 


m  t  x* 

for  any  a  4  H  2  (D  . 


Proof  Using  Green's  formula  in  the  problem  defining  RG,  we 
have 


<iv<H 


)dx  + 


f  i*i: 


(3.22) 


From  (3.22)  it  follows 


(3.23) 

H*ll0,r  <  li«ll0.r 

(3.24) 

•*11*11!  <  II  all 

7  '  r 

Using  the  standard  technique  to  raise  up  the  regularity 
property  of  the  solution  of  partial  differential  equation,  we 
obtain  (3.19)  and  (3.20).  Rewriting  (3.5)  with  an  aid  of  TQ  and 
Re  ,  we  have  for  any  a6Hm+1^2(r) 


(3.25) 


II  R£  (a )  -  a  ||  , 

m- £  ,  r 


✓F  |lTQRC(a)  H  , 

U  m-f,r 


-  0(/e)  ||  RE  (a)  ||  . 

m  +  -j  ,  r 

-  0(/e)  II  a  II  (by  3.19). 

m  + j  ,  r 

Here  we  have  used  the  boundedness  of  Tg  from  Hm+*,'2(r)  to 
(T) . 

Lemma  3.4  Let  e  be  small  enough  and  m£0.  Then 
(3.26)  ||  Se  (a)  +  —  j]  ,  -0(/e))|a||  ,  , 


Proof 


_____  We  prove  this  lemma  in  two  cases.  In  the  first  case, 

2  2 

we  prove  the  special  case;  ■  R+  (flQ«R_)  by  using  the  fourier 
transformation.  We  give  the  plan  of  the  proof  in  the  general 
geometry  in  the  second  case. 


1*  Let  r 


<xx  ,  €)  -  f 

»  «0 


2fiCx 


2  4>e  (Xj^  ,  x2)dx2 


aC5)  -  f  e  15X2  a(x2)dx2,  (a  6  Hn,+3/2(D)  . 

*  —00 


£  '  r 

Here  »J»  (x^,x2)  is  the  solution  of  (3.7)-(3.9).  Then  satisfies 


(3.27) 


i-t-  +  |(l  +  4ir2*|Cj2e)1j/e  «  0  in  R 2  , 
3xf  e  + 


(3.28) 


Xj-0 


-  a  . 


Solving  (3.27)  and  (3.28),  we  have 


(3.29)  <|/e  -  a»exp{-  -(l  +  4Tr2-|C|2*e)7  x,}  . 


from  which 


(3.30) 


fL. 


*i-° 


1 

Se  (a)  -  -  -(1  +  4tt2*  |C|2*e)7  •  a  . 

/e 


We  compute 
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(3.31) 


(3.32) 


S“  (a)  +-§  =  --• - - - - - 

/e  /e  9  o  T 

l+(l+4ir^|cre)^ 


f*°  5  o  n»  -  i  /\  *  2 

<l  +  4ir2-  |c!2)  2  •  Se(e)  +  |  dS 

J  -oo  /e 


<  o(e)  f  |s|2 

*  «vOD  _ 


2-ui2>m+J 


(1+  <1+4tt2*  |  C  1 2e)^]  2 


<0(e)[  |a|2(l  + 4tt2«  U|2)m  +  T  dC 

J  .*00 


which  implies 


(3.33) 


Se  (a)  +  - 


/e  m-i  ,  r 


0(/e)  ||  a  ||  3 

m  +  2  ,  r 


By  density  argument,  we  have 


(3.34) 


Se (a)  +  - 


/e  ||m  -  2  *  r 


0(/e)  ||  a  ||  x 

m  +  ■j  »  r 


m  +  71 


for  va  6  H  ‘(r 


2#  Let  us  now  deal  with  the  general  case.  The  domain  is 
a  regular  simply  connected  domain;  then  there  exists  a  (fixed) 
regular  conformal  mapping  w  »  f  (z)  *  u^  +  iuj  (z  *  x^  +  iXj)  which 


2 

maps  into  B+.  As  a  matter  of  fact,  r  is  mapped  into  the 

£  £  *1 

u2-axis  of  w-plane.  Then  the  transformed  sulution  ¥  »  (f  (w) ) 
satisfies 

12 

H»e  -  0  in  rJ  , 

(3.36)  |^-|  =  A(u2)  =  ^  atf^tw)). 

1  ux-0 

By  means  of  the  iterative  method  proposed  in  the  theory 
of  singular  perturbation  (see  Lions  (1973)),  ipe  is  asymptotically 
developed  in  the  following  ways 

(3.37)  ye  «  I-L  +  y°  +  /e*1  +  ...  +  (/F)n  4*n  +  0n  . 

/c 

Using  (3.37),  we  conclude  (3.26).  ( 


3.4  Let  us  note  that  another  description  of  the  problem  (2.5) 
and  (2.6)  is  following: 

c 

(3.38)  -Aifr0  +  X0<j/Q  *  f  in  nQ  , 

(3.39)  *  0  in  (1^  , 

(3.40)  on  r  , 


(3.41) 


3*0  3*2 

■air  ■  ITT 


on  r  , 
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Then  this  problem  is  transformed  into  the  transmission  equation; 
Find  ae6  Hra+1/2(D  such  that 

(3.43)  Tf (aE)  *  Se(ae). 

Of  course,  the  solution  aE  of  (3.43)  is  equal  to 

3.5  Using  Lemma  3.1,  (3.43)  is' rewritten  in  the  following  way: 


w  m  +  -sr 

(3.44)  Tq (a£-b)  -  Se(ae)  =  -  Tf(b),  for  vb  £  H  (D. 

Recalling  (2.8), 


1 

T 

(3.45)  aE  -  ipe  |  r  -*•  0  strongly  in  H  (D  as  e 0  . 

Therefore  we  choose  b  =  0  in  (3.44).  By  (3.26)  in  Lemma  3.4,  we 
have 


(3.46)  Tn(aE)  +  —  -  SE(aE)  =  -  T-(0). 

0  /e  1  f 

or 

/e«T0(aE)  +  aE  *  /eS^(ae)  -  /£Tf(0) 


where 


e 

S, (aE)  =  SE(aE)  +  —  . 


1 1 


(3.47) 


By  the  definition  of  R  , 


(3.48)  ae  -  /eResJ(ae)  -  /eReTf(0). 


Lemma  3.5  Let  e  be  small  enough.  Then  /eReS^  is  a  contrac¬ 
tion  mapping  from  Hm+*^(r)  onto  itself. 

Proof  Recalling  the  definitions  of  R£  and  S^,  we  see  that 
/cReS£  maps  Hm+1/2(D  onto  Hm+1//2  (T) .  We  compute 


(3.49) 


✓?||ReS£(ae)  ||  ,  <  /e*0 (  ±)||SG(ae)  ||  . 

s  +  »  »  r  /l  m  -  ,  r 


1  ,  n„e  #_£< 


(by  3.20) 

<  0(/e)  ||  aG  ||  , 

i+jj  (by  3.26) 


Thus  (3.49)  completes  this  lemma. 


Using  this  lemma,  we  can  solve  (3.43): 


(3.50) 


(I-/e  ReS^)"1(-/e  ReTf  (0) ) 


B  +  i 

in  H  *(D 


By  (3.21)  in  Lemma  3.3, 


(3.51) 


(I  -  /F  ResJ)"1(-/e  Tf  (0)  +  O (e ) ) 


m  +  § 

in  H  *(r>. 


which  implies 


(3.52) 


in  H 


a  “  -/F«Tj(0)  +  0(e) 


*  ( n . 


Here  we  note  that  Tf(0)  should  satisfy 

0 


Tf<0)  =  -  1TT 


m  +  J 

6  H  iT)  . 


Therefore,  we  have  to  assume  f  6  from  which  follows 


ij;0  e  Hm+3(ft0)  . 


4.  Proof  of  Theorem  3 


4.1  Assume  f €  H  (fig)  (k^5). 
Then 


(4.1) 


<l>°  6  hJ((20)  n  Hk+2(«0)  , 


(4.2) 


e  hK  +  2(r)  n  ck"1,5 


(T) 


(0  <  5  <  1) 


(see  Adams  (1975)).  Using  (2.19)  in  Theorem  3  and  (3.2), 


0(/e)  in  Hk  ?(D  0  Ck_1'5 


(4.3) 


r 


(H 


Using  (4.3)  and  the  maximum  principle,  we  obtain 


(4.4) 


II  *e  llc (ff  )  i  <*>  3)  • 


We  compute 


From  the  definition  of  S  ,  we  have 


k-i 

se(i*»e|r)  e  h  2(d  n  ck“2,6 


By  (3.26)  and  (2.19), 


(4.7) 


^ILt-i.r4*11**  "k  +  i,  r  “  0< 


Combining  (4.6)  and  (4.7), 


(4.8) 


||V4»£  II  i  <0(1) 
k- J,  r 


Similarly,  we  have 


(4.9) 


llVt*  II  i 

k  -  z  ,  3(2 


34i  I 

'k  -  y  ,  3(2 


Since  (4.4),  6H  *(3(2)  H  Ck_2'6  (3(2)  and  <J> 

_  _  9n  3(2 

Put  Vy  81  t  .  Then  H'6  satisfies 


(4.10) 


-eA'T  +  4>  =  0 


in  . 

From  the  maximum  principle  together  with  (4.8)  and  (4.9)  it 
follows 


(4.11)  I|V*€  Ilc(5  }  <  0(1). 


Here  we  have  to  assume  k  >4  to  obtain  good  regularity  of  tj> 
Repeating  similar  argument ,  we  have 


(4.12) 


a2i|/e 


c^) 


<  o(  -) 

/e 


(k  >  5) 


5 .  Proof  of  Theorem  4 

Here  we  restrict  our  proof  to  the  case  XQ  ■  0  without  loss 
of  generality. 


5.1  We  will  introduce  the  following  scalar  product  equivalent  to 
the  one  induced  by  H^(fl): 


(5.1) 


(u  ,  v) 


[  7uW  dx  +  ^  l  oJXvUv)..  h2 

Jn  e  (i,j)en  3 


9 


and  the  associated  norm 


for  any  u,v 


*  hJ(Q) 


(5.2)  ||  u  ||  =  /(u  7  u) 


Let  Hoh(n)={wh6Hj(a) |whl3n  =  0}.  Then,  since  HQh(Q)  is  a 
closed  subspace  in  there  exists  one  and  only  one  element 

eh  e  Hoh(fl)  which  minimize  ||ijje  -  4>^  ”  wh||  among  e  H^h(J2). 

We  denote 


(5.3)  vh  -  *e  -  ♦*  -  eh  . 

Obviously,  satisfies 

(5.4)  (vh  ,  wh)  -  0  , 

and  there  holds 

(s.s)  ||  ♦«  -  ♦«  ||*  -  II  ehH*  +  ||  vh  ||*  . 

Therefore  we  may  evaluate  both  ||  8h  f|  and  |j  vh  |J .  It  follows 
from  (5.4)  that  for  some  C  (see  Ciarlet  (1978)) 

(5.5)  II  vh  II  <  II  *E  -  Ph  II  <  ell  ♦*  ll2,2,n-h 

where  •  J  i|i?  (x)  *  Thus  the  remaining  part  of  this  section 

i  t  j 

will  be  devoted  to  the  estimate  of  |(  8^  || . 


5.2  Using  (5.3),  (5.4)  and  the  relation 


II  «H  H 


VSHitO) 


I <ah  .  v)  | 

II  V  II 


i  <9h :  v  i 

II  wh  || 


(5.7) 


<  Sup 
VvtH^(fl) 


9 


we  have 


I  M>e  -  4>k  ,  w.  )  | 


(5.8)  ||  8.  ||  <  Sup 


h  '  h 


VHoh(n>  11  Wh  11 


Therefore  we  may  compute  (<l>  "  '  wh^  * 

Multiplying  on  both  sides  of  (2.5)  and  integrating  by 
parts,  we  have  (see  (5.1)) 


(5.9) 


(* 


-V  '  L  fKhdx  +  r  J,  ()(h'l'ewh)ij  hi 

»*f  1 1  j 


-I  Us 


dx 


Similarly  from  (2.9)  and  (2.10), 


(5.10) 


(*h'wh>=  l  *  <fwh*^h‘ 
h  h  <i,  j)<?S  »  *5 


oh 


The  subtraction  of  (5.10)  from  (5.9)  gives 


(5.11) 


<*c-*h-v  -  <Lfv-dlt-  (i^an<fVi.j*>2» 


+  {^i?i)^atXh','eWh,i*ih2  ~  M£J,’(’£wh<Sx> 


-  I  +  J 
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Let  Hoh<Q)={wheHh(fi)  !whi  an  =  0}*  Then'  since  HQh(f2)  is  a 
closed  subspace  in  Hj!j(fi),  there  exists  one  and  only  one  element 

E  F  1 

0h  €  Hoh^  which  minimize  ||ip  -  -  wh||  among  wh  e  HQh(ft). 

We  denote 


Obviously,  satisfies 

(5.4)  (vh  ,  wh)  =  0  , 
and  there  holds 

(5.5)  II  -  ^ll2  -  II  «hH2  +  Hvhl|2  . 

Therefore  we  may  evaluate  both  {}  8^  [|  and  ])  vh  jj .  It  follows 
from  (5.4)  that  for  some  C  (see  Ciarlet  (1978)) 


(5.6)  ||  vh  ||  <  ||  <l>6  -  Ph  ||  <  ell  *€  ll2f2,n‘h 

where  p.  **  T  .w^  (x) .  Thus  the  remaining  part  of  this  section 

•  i,j  1,3 

will  be  devoted  to  the  estimate  of  ||  8^  || . 


5.2  Using  (5.3),  (5.4)  and  the  relation 


(5.7) 


ii  «„  ii 


Sup 

v4hJ(0) 


I  (6h  .  v)  | 
II  v  || 


Sup 

VHoh<n> 


I <eh  •  V 1 
II  «K  II 


we  have 


I  (<l»e  -  #  w.  ) 

||  8.  ||  <  Sup  -  h  h 


(5.8) 


VHoh(fi) 


II  wK  || 


Therefore  we  may  compute  '  wh^  * 

Multiplying  wh  on  both  sides  of  (2.5)  and  integrating  by 
parts,  we  have  (see  (5.1)) 


(5.9) 


'"h>  =  L  fwhdx  +  F  A  (*h**ViJ  h: 

“n 


-Hn, 


i|»  w^  dx 


Similarly  from  (2.9)  and  (2.10), 


(5.10) 


(♦h  ,  w)  =  l  (fw.),.^ 

h  h  ti.i)«aoh  h  13 


The  subtraction  of  (5.10)  from  (5.9)  gives 


(5.11) 


(* 


-  .  w.  )  =  (f  fw.  dx  -  l  .  (fw.).  -i  h2^ 

h  h  Jan  h  h  x-3 


+  <r  I  ^(Xk*ew J,  ,h2  -  if  *c  v.  dx > 


Divide  J  into  two  parts  in  the  following  way: 


(5.12) 


J'(|  iW*£wh,i'ih2'^Whdxl 


‘lh 


+  <|  l  a*xh  wh>  i  ih2  "  If  t£whdx) 

e  i,j*&  n  10  e^DrtQ1  n 


J1  +  J2  * 


5.3  In  order  to  estimate  and  J2»  we  will  prepare  the  follow 
ing  lemmas. 


Lemma  5.1  (Quadrature  formula  on  Q) 


(5.13) 


|  [  fw.  dx  -  l  (fwj.  .  h" 

h  (itj)^  h  i'j 


<  ch2  [  .  If  I 


Wh<  ^  I  ^ 


(*> 


(i*j)*?f*  hii 


A1 


where  is  the  support  of  w  J  (see  (2.12). 


Proof:  See  Ciarlet  (1978)  (Error  estimates  for  the  quadrature 
formula  on  a  regular  mesh) 


Lemma  5.2  (Quadrature  formula  on  D) 


(*)  If  I,  _  „  stands  for  the  norm  in  W2,“(o4j). 


(5.14) 


Proof 


|  f  w.dx  -  l  (f  w.  )  .  .  h 

Jd  h  h  13 


<  3/7h  l  |Vf|  *|whiJh 
(i,  '°ij  ftl3 


Since  w.  (x)  =E  w. .w1] (x) 
n  ij  13 


(5.i5,  'Lfwh^-  (irj)ti)«Vijh2i 


l  W(.  (<f  wlj  -  fh2)dx| 
(i.jteO  13  > 


<  i  J“i1ll[  <£»ij-fhV 


There  holds  on  o^: 


<5.16,  I  f  (X)  -  f.J  <  |Vf|  -h/2 

ij  ,0ij 

by  means  of  the  mean  value  theorem.  We  have 

<5.17,  f  wij  dx  -  h2 

and 

2 

(5.18)  area  of  ck  ^  =  3h  . 


Summing  up  (5.15)-(5.18) , 


(5.19)  (f  w*J  -  f  li)dx  <  |7f  |  *h/2*3h  . 

iff .  .  "* 


Substituting  (5.19)  into  (5.15),  we  conclude  (5.14). 


Lemma  5 . 3  (Approximation  of  an  integral  on  a  ship  by  a  boundary 
integral).  With  the  notation  of  the  Figure  4, 


(5.20) 


if  fax  -  f  f  ax| 

JD(\ax  Jq  *i  *2 

<  [  |Vf|  ,  U.  +  l,)2(s)ds 

h  «,rA(si  1  2 


Proof  By  the  mean  values  theorem  for  integrals,  we  have 


(5.21) 


(5.22) 


f  £dx  «  f  <( 

Jr  JrA(s) 


f(x)dr  )dr 


(f  (S)  +Y(s))£1(s)dr 


|y<*>  1  <  |Vf|  •  U.  +  Ij)  (s) 

•,rA(s)  1  2 


By  the  same  argument. 


(5.23) 


if  gar  -  f  7  fy-  a*  I  <  f  |yg|  -a,  + 12)  (s)ds 

Jr  J d  *r*2  Jr  »,rA(s)  1  2 


where  (£^  +  ij)  (x)  is  constant  on  rA  (s)  and  equals  to  (s^  +  Z2)  (s) 


Hence  we  have  (5.20) 


Lemma  5.4  (Discrete  Poincare  inequality) 


(5.24) 


l  wiih2  £  Cf  IVwJ2dx  ' 


where  C  is  the  universal  constant. 


Proof 


(5.25) 


2  „  .  r  ,Wij“WJl-l' 


wij  * 


Prom  (5.25)  , 


(5.26) 


N  a  ^  ~  N  4> .  „  .  .  2  n 

l  "?ih2  <  l  h  (N  l  (—3-^^  )  h2} 

i,j-0  13  ij  4=1 


{ t  )V)M2h2. 

ij  h 


Therefore,  we  have 


(5.27) 


v..2  u2  .  „  **-lk  2  .  ,*ik  ♦ik-l 

4wiih  Ml - H - }  +  ( - h - 

ii  13  Ik  n  n 


Lemma  5.5  (Discrete  trace  theorem) 


(5.28) 


I  *wiih  i  c[  IVwh|2  dx  , 
.,j)f£>  13  'n  h 


Recall  that  (I  a)  < N Z  a 


Let  us  note  that  the  number  of  points  in  D  is  of  order 


Proof 


By  (5.25) 


(5.29) 


T  wt.h  <  [  hN  £  h  ( 

(i,j)*3  “  (i# j)^D 


?  .2  ,w*rw*-ii  2 


<  Cl  | V w.  J 
Jjj  n 


'dx 


Hence  we  conclude  (5.28). 


(i)  Estimation  of  I;  By  Lemma  5.1/ 


(5.30) 


1 

1*1  ll2.-.«ijh2)  X  (J-Wijh2>T 


(ii)  Estimation  of  J^:  By  Lemma  5.1, 


1  1 
2i?/T..2  v.2%7 


(5.31)  Uxl  <  C  11  *e"2,-,oi:jh  »  <Kjh  > 


( iii )  Estimation  of  :  By  Lemmas  5.2  and  5.3, 


(5.32)  |jjl  <  if  *ewh«x  -  f  ♦*'»h  T7nrdx| 

4  Jnnn.  ;D  12 


f  ii;e  w.  .  dx  -  T  „  (Xw  'J,€  W»J  i  ^h2 1 


» 


<  Jrl»«*'VI  ri«*l**2r 


e  *1 


+  ch  uWw«|i™  /h> 


Here  we  see  that  the  definition  of  plays  an  essential  role  in 
this  estimate. 


(iv)  Estimation  of  J:  Therefore 


(5.33) 


|j|  <  %-it  l  ,, 

i.jta, 


1 1  II  2  -  o  1,2)2  (  ?  ?«Wiih2)I 

,0ij  13 


+  (  £  J*<*e  nr’i!  o  h)i  *  I  Xs 

i ,1*0  V*2  ,0ij  i.jtii  13 


?,h>* 


(f  lv*El2  J-<-iH-i>4dr>T-(f  |wh|2  an 
Jr  «,rA  n  Jr  »,r 


<frl*el!(ri-«^)4-)n5-(|rlvwhl2>rlan^ 


By  lemma  5.4  and  5.5,  the  sums  of  in  (5.33)  are  replaced 


by  (/Q|7wh|2dx)1/2. 


Now  iwhL#r.  .  =wkJl,  (kh,*h)  €oij  and  |Vwh|^r  <  2 Iw^^/h 
Therefore,  for  some  curve  ?  near  r 


I a. 


ivw.|  i  2  I  wh «  iul  «  “T  6  «,wn  « 

Jr  h  •,rA  h2  Jr  h  «,rA  h2  (i,j)e?  13 


|  I7W, 


,  C3 

V>"h' 


dx  (by  Lemma  5.S) 


Thus,  summing  up  the  above  results,  we  have 


(5.35)  ||  eh  ||  <  (|I|  +  |j|)/(f  |7wh|2dx) 


1 

7 


<Ch2<  l  .  ||f  ||?  h2) 

<i,j>en  2#  'a<< 


id 


+  %-C(  l  * 

(io)^lh 


h2)7 


e  *1  2 


+  (  I  Jvor  nr-Jl  h) 


i 

,7 


ID 


<f  l**el2  A(^)4dr, 
Jr  »,rA 


l 


♦  %  (f  1*1 


e  ,2  J..+J, 


r1  /e'»,r 


rx,2  4  7 

(-JT~>  dr)2] 


a^TT 


5.5  We  are  now  in  the  final  step: 


(5.36) 


il  r  -  II  =  I!  vh  II  +  II  6h  II 


^  II  *£  II  2  r  2  ,fj 


■h  +  Ch2(  l  ||f||*„0 


+  D(^£). 


(5.37) 


D(^e) 


(  I  *  11/ II2 

(i#  j  )e« 


‘lh 


h2)2 


+  ((Jj)t2) 


e  *1 

v(*  ir-£-) 


*1+Jl2 


h) 


+  <[  |V«l»e|  (-C--  )4  dDJ 

Jr  ®,r  n 


i 


By  using  the  regularity  results  in  Theorem  3, 


(5.38) 


*e  -  *h 


<  c(  -  +  ■£— ) 

l,fl  /e  e/e 


6.  The  Algorithm  for  Solving  Free  Boundary  Problems 

6*1  For  the  sake  of  simplicity,  we  have  tried  our  method  on  a 
simple  geometry  where  the  upper  boundary  of  an  open  rectangle 


is  free  (see  Fig.  5) .  Consider  the  following  problem;  find  $  and 
Y  such  that 


(6.1) 

-  f 

in  fly 

(6.2) 

*  -  0  , 

on  y 

(6.3) 

♦  (0  ,  x2) 

=  4>(1  .  x2) 

(x2  >  0) 

(6.4) 

<Mxx  ,  0) 

=  0 

(0  <  xx  <  1) 

We  assume  that  the  free  boundary  y  is  represented  by  the  equation 
n«n(x.)  (0  <  x^  <  1)  and  satisfy  Q  <  n(x1)  <1  (0  <  <  1} , 


Then  the  discrete  problem  corresponding  to  the  penalized 
state  equation  is  given  by: 


(6.5) 

+  E(xh*C|i,j 

in 

(6.6) 

{0  <  j  <N) 

(6.7) 

<♦‘>1.0  * 

(1  <  i  <N-1)  . 

As  the  boundary  y  (y  ■  n  (x) )  defining  xh»  we  may  choose  the 
linear  interpolation  y^n^Cx)  (0  <  x  <  1)  of  n  *=n(ih)  (O^ij^N). 
We  take  x^  to  be  the  same  as  (2.16)* 

6.2  Define 

(6.8)  g„<x)  -  -i  8‘lx  .  nh(x>) 

(6.9)  g1  -  g(ih  ,  n1) . 


Then  we  introduce  the  discrete  optimization  problem 


(6.10) 


min  (Eh(nh)  *  h J  )gK(ih)  -g1!2)  . 


3E. 

6.3  Owing  to  the  practical  importance  of  »  we  shall  discuss 

3Hh 

its  computation.  For  the  sake  of  simplicity,  we  take  into 
account  only  the  variation  of  y-direction  with  respect  to  differ¬ 
entiation  &t  E.  .  Define 
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(6.11) 


(n  *  nh(ih)) 


n  ■  (ni) 


(6.12) 

n' 

=  (n1  6  5nk*«s 

'Vi 

i. 

Let  6nk  be 

the 

small  deviation 

we  have 

3E. 

N-l 

(6.13) 

n 

V 

“  2h  l  (gh (ih) 

3nx 

i=0  h 

JL  3Vih> 


3n 


Now 

(6.14) 

gh(ih)  =  ±  *h(ih  '  T’h(ih) 

where  6*  is 

defined  by 

(6.15) 

n1  -  hte1^1  -  1)  +  (1  -  e1)  JL] 

Therefore 

(6.16) 

!fh(ih,  *  )ei  + 

3n  /e  3n  i,J  -1 

i#  J 


♦  ^jj($C  i  “  ♦*  £>]  , 

an*  i,J  -l  i,jx 


but  from  (6.15) 

(6.17) 


*  „  v.  36* 

"I? 


To  evaluate  «*» !  .  —  34>f we  differentiate  (6.5)  with  respect  to 


(6.16) 


(Vij  +  F(xh*#)ij 


3n*  h  13  fc  3nk  13 


(plus  some  boundary  conditions  as  in  (6.6)-(6.7) ) 


Let  {p:\}  be  the  solution  of 


(6.i9)  -<v€,ij  +  r(xhp6)ij 

-  (gfc(lh)  -g1)  (0iS  .  +(l-6i)S  J 

h  J.J1-!  j.J1 


p0,j  *  p5.j 

1  <  j  <  N-l 

Pi,«  "  PL.N  -  0 

0  <  1  <  N 

and  multiply  (6.19)  by  and  stun  over  all  i,j;  this  brings 


(6.20) 


l  (ghdh)  -  g1)  (eV  .  +(i-e1)$*  .) 

i  n  i,J  -1  irJx 


ij 


3n 


(J>£)  .  .) 


Therefore 


3E.  N-l  .  ,  i 

(6.21)  — £  -  -2h  l  (g  (ih)  -g1)^ 

3n  i*o  n  3nK 


2h 


I  Ipf.(A<XK«..+r( 


3x 


h  j.c 


n  ij  1 «  »nk,Ah  1 


3n 


3n 


13 


2  / _  /i.L  \  _k %  i .  £ 


i<g„(kh>  -3“)  (4‘  -*E  „>  =  M* 


(6.11) 


(n  «  nh(ih)) 


n  -  (n  } 

<\# 


(6.12) 

n*  -  (n1  6  <$nk*<5 .  ,  )  . 

*  i,k 

Let  6nk  be 

the  small  deviation  of  the  k  th  component  of  n. 

Then 

we  have 

(6.13) 

3E.  N-l  .  3g.  (ih) 

l.  “  2h  l  (g.  (ih)  -g  )  - j- —  . 

3n  i=o  n  3n 

Now 

(6.14) 

g.(ih)  -  i  $5  (ih  ,  nK(ih)  -  -  [  4>e  .e1  +  4e  ,  (1 

h  /l  h  h  /e  i,JL  i.J1 

-  0*)] 

where  0i  is 

defined  by 

(6.15) 

n1  -  hie1^1  -  1)  +  (1  -  eA)  J1! 

Therefore 

(6.16) 

!fh(ih,  «  3  )0i  +  -i_ (4>e  i)(l-0i) 

3n  /e  3n  i,J  -1  3n  i,J 

.  30  ..6  ,6  »  i 

+  i  ”  ♦  i >1  » 

3n  i,J  -1  i,J 

but  from  (6 

.15) 

(6.17) 

*  -  K3®1 

®  i  v  *  r  . 

To  evaluate  » 3^i/3nk  we  differentiate  (6.5)  with  respect  to 


(6.18) 


(Vij  + 


8  /  1 ,8xh  Ae» 

^lc  b  i3  e  3nk  13 


(plus  some  boundary  conditions  as  in  (6.6)— (6.7) ) 


Let  tp^j)  be  the  solution  of 


(6.19) 


-<V€,ij  +  r()thpC>ij 


(g.(ih)  -  g1)  (0i6  *  +(l-0i)6  <> 


j,j -i 


3»J 


p0,j  -  PN,j 


Pi,0  "  PL,N 


1  <  j  <  N-l 


0  <  1  <  N 


and  multiply  (6.19)  by  <J>!^  and  sura  over  all  i,j;  this  brings 


(6.20) 


Therefore 


(6.21) 


£(g  (ih)  -g1)  (eV  +  (1  -  O1)#'  .) 

i  n  i,J  -1  i^J1 


’  2  pi j  (“1?  (xhf }  ii  +  F(“T  )  ii > 
ij  13  3n*  n  13  e  3nk  13 


BE.  n-i  .  i 

— $  »  -2h  l  (gh(ih)  -g1)^- 
3n  i=0  n  3n* 


.  |(gh(kh)  -gk)  (*e  -  <pe  k)  i  Mk 

/l  h  k,J-l  k,JK  h 


*1 


-4 

•I 
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(6.22) 


Take  6n  «  -  A  MjJ  (0  <  A  <<  1) 


Then  we  have 


(6.23)  Eh(Jj’}  “  Eh<J)  "  “  <  0. 


6.4  Therefore  let  us  consider  the  following  algorithm  for 
solving  (6.10): 

step  1  Choose  A  >  0  small,  choose  the  initial  guess  of  the 
free  boundary 

% 

step  2  Determine  xk  by 

step  3  Compute  by  solving  (6.5)- (6.7). 
step  4  Compute  p^  by  solving  (6.14)- (6.16) . 
step  5  Set 

6nk  -  -  AMk  (0  <  k  <  N-l) . 

and 

-  {ni  +  6nk  •  6ki) 

replace  nh  by  and  go  back  to  step  2. 

6.5  As  an  example,  we  shall  deal  with  the  problem  (6.1)-(6.4). 
Which  has  a  solution  (y  ■  n  (x)  *  a  +  b  sin  2irx,  <Mx,y)  -  cy  (n  (x)  -  y) } 
Naturally,  the  data  is  defined 


(6.24) 


C-  (2  -  y  2_£  )  =  C  (2  +  4tt2  y  sin2irx) 
dx 


(6.25) 


-  c  n 


l  l 

t1  +  (gj)  2)^  »  -  C  sin  2ttx  (1  +  4n2  sin2  2*)^ 


6.5  By  using  the  algorithm  mentioned  above,  we  have  obtained 
several  numerical  results,  which  will  be  reported  in  following 


'-Si 


papers. 
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